﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using Meta.Numerics.Statistics.Distributions;

namespace ExcelAddIn1
{
    /**
     * Implementation of the Wilcoxon rank-sum test.
     */
    class Test_WilcoxonSumRank : Test
    {
        private class RankItem : IComparable
        {
            public double value;
            public String dataset;
            public double rank;
        
            public RankItem(double value, String dataset)
            {
                this.value = value;
                this.dataset = dataset;
                this.rank = 0;
            }
        
            public int CompareTo(Object obj)
            {
                RankItem ri = (RankItem)obj;
                if (value < ri.value)
                {
                    return -1;
                }
                else if (value > ri.value)
                {
                    return 1;
                }
                else
                {
                    return 0;
                }
            }
        }

        private String[,] Rc005 = {
            {"5;16", "6;18", "6;21", "7;23", "7;26", "8;28", "8;31", "9;33"},
            {"6;18", "11;25", "12;28", "12;32", "13;35", "14;38", "15;41", "16;44"},
            {"6;21", "12;28", "18;37", "19;41", "20;45", "21;49", "22;53", "24;56"},
            {"7;23", "12;32", "19;41", "26;5", "28;56", "29;61", "31;65", "32;70"},
            {"7;26", "13;35", "20;45", "28;56", "37;68", "39;73", "41;78", "43;83"},
            {"8;28", "14;38", "21;49", "29;61", "39;73", "49;87", "51;93", "54;98"},
            {"8;31", "15;41", "22;53", "31;65", "41;78", "51;93", "63;108", "66;114"},
            {"9;33", "16;44", "24;56", "32;70", "43;83", "54;98", "66;114", "79;131"}
        };

        private String[,] Rc010 = {
            {"6;15", "7;17", "7;20", "8;22", "9;24", "9;27", "10;29", "11;31"},
            {"7;17", "12;24", "13;27", "14;30", "15;33", "16;36", "17;39", "18;42"},
            {"7;20", "13;27", "19;36", "20;40", "22;43", "24;46", "25;50", "26;54"},
            {"8;22", "14;30", "20;40", "28;50", "30;54", "32;58", "33;63", "35;67"},
            {"9;24", "15;33", "22;43", "30;54", "39;66", "41;71", "43;76", "46;80"},
            {"9;27", "16;36", "24;46", "32;58", "41;71", "52;84", "54;90", "57;95"},
            {"10;29", "17;39", "25;50", "33;63", "43;76", "54;90", "66;105", "69;111"},
            {"11;31", "18;42", "26;54", "35;67", "46;80", "57;95", "69;111", "83;127"}
        };

        public Test_WilcoxonSumRank(DataContainer data)
        {
            this.data = data;
            res = new List<String>();
        }

        public override bool EnableAlphaField()
        {
            return false;
        }

        public override string GetInfo()
        {
            String s = "Performs a Wilcoxon rank-sum test to check for significant differences between two datasets. ";
            s += "The test does not require normally distributed data, but the datasets must be independent. This test is also ";
            s += "known as the Wilcoxon-Mann-Whitney test.";
            return s;
        }

        public override void RunTest()
        {
            res.Add("Wilcoxon rank-sum test");

            //Error check
            if (data.GetNoSets() != 2)
            {
                res.Add("Two datasets is required!");
                return;
            }
            if (alpha != 0.05)
            {
                if (alpha != 0.1)
                {
                    res.Add("The test equires Alpha=0.05 or Alpha=0.10");
                    return;
                }
            }

            //Get data sets
            DataSet d1 = data.GetDataSet(0);
            DataSet d2 = data.GetDataSet(1);
            int totalN = d1.GetN() + d2.GetN();

            //Create the rankings table
            List<RankItem> rankTab = new List<RankItem>(totalN);
            //Fill table with both datasets
            for (int i = 0; i < d1.GetN(); i++)
            {
                rankTab.Add(new RankItem(d1.GetValue(i), d1.GetName()));
            }
            for (int i = 0; i < d2.GetN(); i++)
            {
                rankTab.Add(new RankItem(d2.GetValue(i), d2.GetName()));
            }
            //Sort the table
            rankTab.Sort();
            //Fill in ranks
            for (int i = 0; i < totalN; i++)
            {
                rankTab[i].rank = i + 1;
            }
            //Fix rank duplicates
            double prev = rankTab[0].value;
            double sumRank = rankTab[0].rank;
            int cnt = 1;
            for (int i = 1; i < totalN; i++)
            {
                double iVal = rankTab[i].value;
                if (iVal == prev)
                {
                    //Equal values. Find a mean rank.
                    sumRank += rankTab[i].rank;
                    cnt++;
                }
                else
                {
                    if (cnt > 1)
                    {
                        //Duplicates found. Update rank to a mean rank.
                        double meanR = sumRank / (double)(cnt);
                        for (int d = i - 1; d >= i - cnt; d--)
                        {
                            rankTab[d].rank = meanR;
                        }
                    }
                    cnt = 1;
                    sumRank = rankTab[i].rank;
                }
                prev = iVal;
            }
            //Rankings table done!

            //Calcualte mean rank for the dataset
            //with fewest samples.
            DataSet rankData = d1;
            if (d2.GetN() < d1.GetN())
            {
                rankData = d2;
            }

            double R = 0.0;
            for (int i = 0; i < totalN; i++)
            {
                RankItem ri = rankTab[i];
                if (ri.dataset.Equals(rankData.GetName()))
                {
                    R += ri.rank;
                }
            }

            int i1 = d1.GetN() - 3;
            int i2 = d2.GetN() - 3;

            //For low sample size, use the Wilcoxon table
            if (d1.GetN() <= 10 && d2.GetN() <= 10)
            {
                //Find critical R
                String[] tokens = null;
                if (alpha == 0.05)
                {
                    tokens = Rc005[i1,i2].Split(';');
                }
                else
                {
                    tokens = Rc010[i1,i2].Split(';');
                }
                int RcLow = Convert.ToInt32(tokens[0]);
                int RcHigh = Convert.ToInt32(tokens[1]);

                //Results
                res.Add(";" + d1.GetName() + ";" + d2.GetName() + ";-");
                res.Add("N;" + d1.GetN() + ";" + d2.GetN() + ";-");
                res.Add("Mean;" + d1.GetMean().ToString("F2") + ";" + d2.GetMean().ToString("F2"));
                res.Add("StDev;" + d1.GetStDev().ToString("F2") + ";" + d2.GetStDev().ToString("F2"));
                res.Add("α;" + alpha + ";;-");
                res.Add("R-score;" + R.ToString("F2"));
                res.Add("R-crit;" + RcLow.ToString("F2") + ";" + RcHigh.ToString("F2"));
                if (R < RcLow || R > RcHigh)
                {
                    res.Add("Result;Significant difference");
                }
                else
                {
                    res.Add("Result;No difference");
                }
                res.Add(";;;-");
            }
            else
            {
                //For large sample size, approximate using normal distribution.
                DataSet l1 = d1;
                DataSet l2 = d2;
                if (d1.GetN() > d2.GetN())
                {
                    l1 = d2;
                    l2 = d1;
                }
                double n1 = (double)l1.GetN();
                double n2 = (double)l2.GetN();

                double mean = n1 * (n1 + n2 + 1.0) / 2.0;
                double stdev = Math.Sqrt(n1 * n2 * (n1 + n2 + 1.0) / (Math.Max(n1, n2)));

                double Z = (R - mean) / stdev;

                NormalDistribution nDist = new NormalDistribution();
                double ZcLow = nDist.InverseLeftProbability(alpha / 2.0);
                double ZcHigh = nDist.InverseRightProbability(alpha / 2.0);

                //Results
                res.Add(";" + d1.GetName() + ";" + d2.GetName() + ";-");
                res.Add("N;" + d1.GetN() + ";" + d2.GetN() + ";-");
                res.Add("Mean;" + d1.GetMean().ToString("F2") + ";" + d2.GetMean().ToString("F2"));
                res.Add("StDev;" + d1.GetStDev().ToString("F2") + ";" + d2.GetStDev().ToString("F2"));
                res.Add("Alpha;" + alpha + ";;-");
                res.Add("Z-score;" + Z.ToString("F2"));
                res.Add("Z-crit;" + ZcHigh.ToString("F2"));
                if (Z < ZcLow || Z > ZcHigh)
                {
                    res.Add("Result;Significant difference");
                }
                else
                {
                    res.Add("Result;No difference");
                }
                res.Add(";;;-");
            }
        }
    }
}
